Homogeneous isotropic turbulence in dilute polymers: scale by scale budget 
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ABSTRACT 

The turbulent energy cascade in dilute polymers so- 
lution is addressed here by considering a direct numeri- 
cal simulation of homogeneous isotropic turbulence of a 
FENE-P fluid in a triply periodic box. On the basis of 
the DNS data, a scale by scale analysis is provided by 
using the proper extension to visco-elastic fluids of the 
Karman-Howarth equation for the velocity. For the mi- 
crostructure, an equation, analogous to the Yaglom equa- 
tion for scalars, is proposed for the free-energy density 
associated to the elastic behavior of the material. Two 
mechanisms of energy removal from the scale of the forc- 
ing are identified, namely the classical non-linear trans- 
fer term of the standard Navier-Stokes equations and the 
coupling between macroscopic velocity and microstruc- 
ture. The latter, on average, drains kinetic energy to feed 
the dynamics of the microstructure. The cross-over scale 
between the two corresponding energy fluxes is identified, 
with the flux associated with the microstructure domi- 
nating at small separations to become sub-leading above 
the cross-over scale, which is the equivalent of the elastic 
limit scale defined by De Gennes- Tabor on the basis of 
phenomenological assumptions. 



I. INTRODUCTION 

Turbulence in dilute polymers solutions is still an open 
issue, despite the growing interest on the subject. Most 
of the efforts have been addressed to the comprehension 
of practical aspects, such as drag reduction, related to 
the modification of turbulence by long chain polymers in 
wall bounded flows, |l7j ]. However, in most cases, the 
attempts have proven inconclusive, suggesting that, even 
from the point of view of applications, a more fundamen- 
tal approach is required. Under this respect, basic studies 
should focus on simpler flow conditions where turbulence- 
polymers interaction occurs. 

A possible choice is to exploit homogeneity and 
isotropy to simplify the statistical treatment of the data 
and to make use of the available, relatively well estab- 
lished, rheological models for dilute polymers, see Q for 
a comprehensive review. 
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In this context, a recent experimental investigation on 
decaying grid turbulence of dilute polymers and surfac- 
tants [jlQ] shows, for the polymers, a substantial alter- 
ation of the small scales, though still in presence of a sub- 
stantial degree of anisotropy, ascribed to the grid shear 
layers. 

In this and other related experiments, see e.g. j23] and 
[ p"l| , the measurements were mostly aimed at the velocity, 
but also at the pressure field see e.g. |3| for a recent exam- 
ple. Actually, as one can imagine, the experimental anal- 
ysis of turbulent flows of polymeric liquids is particularly 
difficult, and hardly one can proceed beyond the global 
characterization of the flow and the analysis of certain as- 
pects of the macroscopic field, see Jig] , [ p6| as additional 
examples in the more complex configuration of a channel 
flow. This makes the actual mechanism of interaction be- 
tween polymers and turbulence particularly obscure, and 
leads to consider numerical simulations as a viable tool to 
address the problem. For triply-periodic boundary con- 
ditions, in particular, fast and accurate spectral methods 
can easily be developed to achieve the direct numerical 
simulation of a strictly homogeneous isotropic flow, once 
a reasonably accurate and sufficiently simple rheologi- 
cal model has been selected. Among these, the so-called 
FENE-P model 25 is possibly the best compromise be- 
tween accurate representation of polymers dynamics and 
minimal computational complexity, see pTEj ] . In fact, re- 
cent results show that this model is able to reproduce 
the drag reducing behavior of dilute polymers in wall 
bounded flows, see §|, §, §. 

As often in rheology, one can easily get stuck in a dis- 
cussion about the proper model for a specific application. 
We take here the model for granted and try to grasp the 
mechanisms by which it affects the turbulence. To do 
this we take freedom to chose parameter values which 
may not correspond too well to experimental conditions, 
but which may help revealing the underlying physics. 

In summary, our model for the solution considers an 
ensemble of elastic dumbbells attached to each material 
point of the continuum. The dumbbells are stretched by 
the flow due to friction and react through a non-linear 
clastic spring. The probability density function of the 
end-to-end vector of the dumbbells is described through 
a second order tensor, the conformation tensor, assumed 
to represent the covariance matrix of the end-to-end vec- 
tors. It obeys a transport equation which accounts for 
material stretching and for the elastic reaction. The force 
exerted on the continuum by the dumbells contributes an 
extra-stress term in the macroscopic momentum balance. 
On average, the work done by the extra-stress against 
the macroscopic velocity field amounts to a draining of 
energy from the macroscopic field. The corresponding 
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power feeds the dynamics of the microstructure, is par- 
tially accumulated as free-energy of the ensemble and is 
eventually dissipated, see ||. 

In terms of turbulence dynamics, several basic ques- 
tions arise. First of all, one should understand if the 
nature of this kind of viscoelastic turbulence is substan- 
tially similar, or, on the contrary, essentially different 
from standard Newtonian turbulence. According to the 
common view, a turbulent flow is described as a super- 
position of different scales of motion. In homogeneous 
isotropic turbulence of ordinary fluids the interaction of 
the different scales originates a flux of energy from the 
large towards the small scale, with energy injected by 
an external forcing at large scales and dissipated by vis- 
cosity at small scales, the so-called direct, or forward, 
energy cascade (see e.g. the book by Frisch Jl2|]). The 
presence of the microstructure is expected to alter this 
process, since it corresponds to an alternative path the 
energy flux may take to achieve the eventual dissipation. 
If this is the case, one may wonder if the alteration is 
localized at certain spatial scales or if the entire spec- 
trum is affected by the energy draining of the polymers. 
This poses the additional issue of determining a possible 
scale which delimits the range of scales where the pres- 
ence of the polymers is effective from those where they 
are substantially passive, see ||. This scale may be close 
to the classical Kolmogorov dissipation scale, but it may 
even be substantially larger, and it is not even clear if 
the effect of polymers may be reduced to an additional 
dissipation, |lq], or if the role of elastic energy is crucial, 



see also 23 



All these issues can be addressed systematically by an- 
alyzing the data of a DNS of homogeneous isotropic tur- 
bulence once proper diagnostic tools are available. To 
this purpose we present and discuss a scale by scale bud- 
get based on the extension to viscoelastic fluids of the 
classical Karman-Howarth Jl2| equation of ordinary tur- 
bulence, see e.g. [Q, |^] and j|] for similar extensions 
in the context of shear flows. By this equation we can 
address the effect of polymers on the macroscopic kinetic 
field. To consider the fluctuations in the polymers we 
would need to extend the classical Yaglom equation for 
a scalar, see ^f|, to a second order active tensor. In 
a more straightforward way instead we consider, in the 
proper thermodynamic framework, the Yaglom equation 
for the free-energy of the polymers. 



II. THE EVOLUTION EQUATIONS FOR DILUTE 
POLYMER SOLUTIONS 

As discussed in the introduction, the momentum bal- 
ance for a dilute solution of long chain polymers in oth- 
erwise Newtonian, incompressible solvent is described by 
a slightly modified form of the Navier-Stokes equation, 
namely 



~~Dt 



dp 

dxi 



d 2 u, 
dxjdxj 



dxj 



+ fi, 



(1) 



where D/Dt = d/dt + Ukd/dxk is the usual substan- 
tial derivative, the external forcing, p is the pressure 
normalized by the constant density of the solution, v is 
the kinematic viscosity of the solvent and T]y is the ad- 
ditional contribution to the stress tensor - extra-stress - 
due to the polymers chains. The macroscopic velocity 
field Ui is solenoidal and the constitutive relation for the 
extra-stress is 



Tij = v v /r [f{Rkk] Pm, Po)Rij I 'Po 



(2) 



where v p is a constant typically of the order of a small 
but finite fraction of v. The conformation tensor, Rij, 
is a second order tensor taken to characterize the sta- 
tistical behavior of the population of polymers chains 
attached to a given point, and f{Rkk',Pm,Po) is a di- 
mensionless nonlinear spring coefficient specified as / = 
(Pm ~ Po)/(Pm ~ Rkk), with po the equilibrium length 
of the chains and p m their maximum allowed length, as 
standard in the context of the finite extensibility nonlin- 
ear elastic model (FENE, see ||). The dynamics of the 
population of chains is described by an evolution equa- 
tion for the conformation tensor, 



DR 



Dt 



duk 

dx-j 



(3) 



where r is the principal relaxation time of the chains. 
The system (Q, ||, ||) completed with the continuity equa- 
tion, 



dxi 



0. 



(4) 



for the macroscopic velocity forms the FENE-P model 
for dilute polymers solutions, see e.g. ^BJ. 
Phenomenologically, stretching of the chains by the flow, 
described by the terms involving the velocity gradient 
in (|^), is counteracted by an elastic reaction, the term 
in brackets, while the polymers are transported and re- 
oriented by the flow. The reaction of the polymers alters 
the macroscopic force balance, as seen in eq. (|l|), where 
the divergence of the extra-stress enters as an additional 
source or sink of momentum. 

The extra-stress is able to make work against the ve- 
locity as follows from the balance of kinetic energy 



D(u 2 /2) dipSij + E^m , d(T l]Ul ) 



Dt 



dXj 



dXj 



e N + S, (5) 



where = 2vaj is the standard viscous component 
of the stress tensor and ejy = 21^6™ ey is the Newtonian 
component of the dissipation, with e.y the macroscopic 
velocity deformation tensor. With respect to the kinetic 
energy balance of the standard Navier-Stokes equation, 
two additional terms arise here, to describe the inter- 
action of the macroscopic field with the microstructure. 
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The one in divergence form contributes to the spatial re- 
distribution of kinetic energy. The other, 



S = T t 



duj 



(6) 



called hereafter the stress-power, represent energy per 
unit time that the microstructure can drain, or, in prin- 
ciple, even release to the macroscopic kinetic field. The 
presence of this term is physically due to the fact that the 
polymers chains can both dissipate, via friction with the 
solvent, or store energy, due to their elasticity. The cor- 
rect thermodynamical setting refers to the free-energy a 
of the polymers ensemble. Given the expression, see e.g. 



1 Up 

2 T 



(pLax/Po - 1) log 



R, 



-Pa 



-1/3 \og(det R/pD] 



(7) 



it follows directly from the evolution equation (||) that 
the energy balance for the polymers is 



Da 

Dt ~ b 



where 



1 v. 



(8) 



(9) 



a positive definite quantity, is the rate of energy dissipa- 
tion associated with the polymers. Adding eq. (||) and (||) 
gives the balance of the total free-energy of the system, 
E = u 2 /2 + a, as 



DE dJ k 
-prr = fiUi -e T + ^~ 
Dt oxk 



(10) 



where Ji is the spatial flux of total free-energy, whose ex- 
pression follows immediately by inspection of the kinetic 
energy balance. Note that S drops from the global bal- 
ance since it corresponds to a conservative exchange of 
energy between the two constituents, namely the macro- 
scopic kinetic energy and the polymers free-energy. 



III. THE KARMAN-HOWARTH EQUATION 

In a statistically stationary turbulent system, the ex- 
ternal forcing fa provides the energy to maintain the tur- 
bulent fluctuations which otherwise would be damped by 
the dissipation. For the present case both the macro- 
scopic velocity and the microstructure exhibits turbulent 
fluctuations. Characteristic feature of the system is that, 
while the external forcing acts only on the macroscopic 
field, both the macroscopic flow and the microstructure 
contribute to the dissipation of the injected mechanical 
power. This is possible since the extra-stress captures 



part of the available mechanical energy to feed the fluc- 
tuations of the microstructure and the associated dissi- 
pation. In order to try to disentangle these complex pro- 
cesses, the simplest setting is provided by homogeneous 
and isotropic turbulence once suitable balance equation 
are available for the fluctuations of both kinds. 

Purpose of the present section is to present the equa- 
tion appropriate to address the scale by scale budget of 
the macroscopic velocity fluctuations. This is the ex- 
tention to dilute polymer solutions of the well known 
Karman-Howarth equation of homogeneous isotropic tur- 
bulence for Newtonian fluids. As discussed in the previ- 
ous section, where the FENE-P model has been briefly 
recalled, the solenoidal velocity field, Ui, obeys a modified 
form of the Navier-Stokes equation, where the divergence 
of the extra-stress T, denoted here by 



(11) 



appears as an additional source of momentum. 

Starting with the modified Navier-Stokes equation, the 
equation for the correlation tensor Cij —< UiU 1 ^ >, 
where unprimed and primed variables are evaluated at 
Xi and at x' i} respectively, with x\ the displacement of 
xi by the separation vector r^, easily follows from the 
standard procedure, see e.g. Q, The trace of this 
equation reads 



2v 



dCj,i 
dt 
d 2 Ci,j 
dr k dr t 



d 
dr k 



< UiU { U k - UiUiUk >- 



+ < g\ui + gtu'i > + < f[ui + > 



(12) 



where fi is here specified as a random stirring force acting 
on the large scales. To obtain the above equation we have 
used standard properties of homogeneous incompressible 
turbulence, which allow to drop the pressure term and to 
recast spatial derivatives with derivatives with respect to 
the separation, see e.g. Q. 

Equation ( |l2| ) can be re-expressed in terms of velocity 
increments, 5Vi = Ui(x' r ) — Ui(x r ), SV 2 ~ SViSVi, once 
recalled that C^i =< Uim > —1/2 < SV 2 >, and that, 
by standard manipulations, we have 



< 8V 2 SV k >= -2 < Uiu'^ - Utii'tUk > + 

< UiUiu' k > — < u'iU'iUk > ■ 



(13) 



Hence, for a solenoidal velocity field, the equation in 
terms of velocity increments follows as 



d 
dr k 



< SV 2 SV k + 2T* kl W t >= -4 < f lUl > 



2 < SVi5fi > + 2v 



d 2 



dr k dr k 



< SV 2 > 



(14) 
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In deriving eq. ( JL4| ) a statistically stationary field has 
been assumed and the term < f[ui + fiU^ > has been re- 
arranged in the form 2 < /jitj > — < SViSfi > to extract 
explicitly the power injected by the external forcing. Fi- 
nally, the integration by parts 



< g'i^i + giu't >= -^-{T' ki Ui 

UT}- 



d 



where = TL + T kl , leads to eq. @. 

Equation (y_4j) can be integrated in r-space over a ball 
B r of radius r (see e.g. jlj for a similar approach) to yield 

^-2 <f> < 8V 2 8V k + 2 T ki SVi > n k dS r = -4/3 e T r 



27rr 2 dr 
1 



— * < <SU 2 > dS* r 



for the velocity field. Its Fourier transform actually cor- 
responds to the budget of turbulent kinetic energy in the 
different modes of the spectrum. To address the complete 
dynamics of dilute polymers solution we need a similar 
equation concerning the microstructure. Under this re- 
spect, in section 1 we have shown that a free-energy a can 
be associated to the ensemble of polymers attached to a 
given point in space and that, locally, the energy balance 
for the microstructure is provided by equation (||). 

The free-energy a is a measure of the average elonga- 
tion of the polymers as given by the trace of the confor- 
mation tensor Rij . In fact one can derive an equation for 
the scale by scale budget of the trace though we prefer 
to discuss the Yaglom equation for a, whose meaning in 
terms of energetics of the system is more definite. To this 
purpose, let us denote by 



dB r 



2irr 2 



1> 



(18) 



where the global energy balance 
tionary homogeneous field reads 



< SVififi > dV r 



which for a sta- 



(17) 



< fiUi > 

has been accounted for. 

The Karman-Howarth equation ( |l6| ) highlights the 
main difference between Newtonian and viscoelastic tur- 
bulence. For either case, in the range of scales where 
the direct effect of viscosity and the correlation between 
velocity and forcing increments are both negligible, the 
energy flux occurring through scale r equals the total 
dissipation. For Newtonian fluids, all the dissipation is 
provided by the viscosity and the flux through scale r 
only occurs due to the classical nonlinearity associated to 
the advection. For viscoelastic fluids, an additional dis- 
sipative process takes place in the polymers, and, consis- 
tently, the flux through scale r also presents a new com- 
ponent in addition to that of the standard Navier-Stokes 
dynamics. In fact, the additional flux term describes 
the amount of energy intercepted by the microstructure. 
Though the classic component is expected to be associ- 
ated to a forward energy cascade, i.e. the corresponding 
flux enters the ball of radius r, concerning the viscoelas- 
tic contribution we can, at the moment, only conjecture 
that energy intercepted at scale r is passed forward to 
the microstructure at smaller scale, consistently with the 
idea of polymer dissipation as a small scale process. How- 
ever, since, to our knowledge, the process has never been 
investigated before, we have to examine and discuss the 
numerical simulation before drawing any conclusion. 



IV. THE YAGLOM EQUATION FOR THE 
FREE-ENERGY 

The Karman-Howarth equation for dilute polymer so- 
lution may be used to analyze the scale by scale budget 



the instantaneous excess power, i.e. the amount of power 
(15) transfered to the polymers in excess to the dissipation 
they originate. Strictly speaking this is an improper 
naming, since on the average the excess power is zero. 
Nonetheless, it measures the amount of energy stored 
per unit time in the polymers. From (||) we can derive 
the equation for the correlation < aa' >, 



d < aa' > 
di 



d 
dr k 



< u k a a — Ukaa > 
< a'lp + atp' > , 



(19) 



which, considering that 
< atp' + a'lp >= 



2 < aip > - < SaSi/j > , 



and noting that, for a stationary state (< aDa/Dt >- 
0), 



< aip > 

can be recast in the form 
i) 



0. 



dr k 



< Sukda- 2 >= 2 < SaSip > 



(20) 



(21) 



Integration over the sphere of radius r then yields 



1 



^r 2 JdB , 
2 



< Sa 2 5V k > n k dS r 



< SaSip > dV r , 



47rr 2 j Br 

which, by using isotropy reduces to 



< 8a 2 5V\\ > = 2 f 
Jo 



< SaSrp > r 2 dr , 



(22) 



(23) 



with <5Vj| = SVi Tijr the longitudinal velocity increment. 

Equation (^3|) is the appropriate form of the Yaglom 
equation for the polymers. It states the equivalence be- 
tween the flux of free energy through the different spatial 
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scales, as described by < 5a 2 5V\\ > and the integral of the 
correlation between free energy fluctuations, 5a, and fluc- 
tuation of the excess of power, Sip. Classically for scalars 
the Yaglom equation provides a term linear in the sepa- 
ration associated to the average scalar dissipation. In the 
present case, such term is zero as consequence of eq. (20). 



V. GLOBAL FEATURES OF THE FLOW AND 
DNS 

The system of equations (|l],|]j4]), corresponding to 9 
scalar equation, 3 velocity components plus 6 compo- 
nents of the conformation tensor, has been integrated 
using a Fourier x Fourier x Fourier spectral method 
for spatial discretization and a four-stages/third order 
Runge-Kutta method for time advancement with all the 
non-linear terms fully de-aliased by the three-halves rule. 
An exactly solenoidal velocity is achieved by the stan- 
dard projection method. The side of the cubic compu- 
tational box is 2tt and 64 modes are used in each direc- 
tion. A fully developed field of Newtonian turbulence has 
been used as initial condition for the viscoelastic calcu- 
lation and the random forcing is applied to the first shell 
of wave-vectors, with constant amplitude and uniformly 
distributed phases. The random phase extracted at each 
step is used subject to the condition of positive energy 
injection, otherwise it is changed by 7r. 

The nominal Reynolds number of the simulation can 
be constructed from the intensity of the forcing, /o, the 
length of the box, l x , and the zero shear kinematic vis- 
cosity, v T = v + up, as Re = (^/o) (1/2) /^t = 676. The 
parameters concerning the polymers are r\ v = v v jv = .1, 
see e.g. @, p 2 J p 2 = 1000 and De = r(/ // x ) 1 /2 = 
.38. When comparing viscoelastic and Newtonian simu- 
lations, we keep constant the nominal Reynolds number, 
i.e. for the Newtonian case we use the same forcing in- 
tensity and a value for the viscosity equal to the value of 
vt for the viscoelastic fluid. 

After the transient is elapsed, statistically stationary 
conditions are achieved. The general impression on the 
field may be gained by considering figure [l], where instan- 
taneous vortical structures are shown for the Newtonian 
and the viscoelastic simulation. The eduction criterion 
for vortical structures is based on the discriminant of 
the velocity gradient tensor, one of the now classical ap- 
proaches to extract vortical structures from a complex 
field, see e.g. || for details. Algebraically, the idea is 
to identify with vortical structures those regions where a 
couple of eigenvalues of the velocity gradient are complex 
conjugate, implying a locally helicoidal relative motion 
of nearby particles. In principle, other criteria could be 
employed and different field quantities could be visual- 
ized. In any case, the inspection of the instantaneous 
fields always leads to the conclusion that, with respect to 
Newtonian, much less but wider structures are present in 
the viscoelastic case. 




FIG. 1. Vortical structures for the Newtonian, De = 0, 
(top) and viscoelastic case, De = .38, (bottom). 

Table | presents the global parameters of a typical 
viscoelastic simulation, De — .38, together with a cor- 
responding Newtonian case. In the table, the Taylor 
microscale is defined as A = (l5/2u 2 ms /Q) 1 / 2 , with 
n = l/2<( 2 >,( being the vorticity. The average value 
of the free-energy is denoted by a, while the other sym- 
bols have been already defined in the previous sections. 
Concerning the Kolmogorov scale, r\, for the viscoelastic 
case this quantity is not uniquely defined. Actually, for 
Newtonian turbulence one assumes that the small scale 
dynamics is controlled by dissipation and viscosity, orig- 
inating the classical definition of the Kolmogorov length 
as rj = (i^/e) 1 / 4 . In the present case, several other pa- 
rameters related to the more complex rheology may in 
principle affect the small scale dynamics. In the table we 
have conventionally assumed t\t = (^tAt) 1 ^ 4 - 



De 


Rex 




Q 


a £t 




ep A 


Vt 





75 


.806 


10.6 


- .156 


.156 


.68 


.040 


.38 


170 


.9685 


5.1 


.09 .231 


.067 


.164 1.18 


.036 



TABLE I. Global parameters for the Newtonian, De = 0, 
and the viscoelastic, De = .38, simulation. 
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By discussing the data in table Q, the macroscopic ef- 
fects of viscoelasticity can be summarized as follows. For 
fixed forcing amplitude, i.e. for fixed nominal Reynolds 
number in identical boxes and with same total viscosity, 
Vti we observe a substantial increase of the total dis- 
sipation. Clearly, in the steady state, this corresponds 
to an equivalent increase in the power extracted from 
the external forcing. For the viscoelastic case, most of 
the dissipation, about 71%, occurs in the polymers. The 
u rms increases of about 20%, with respect to the New- 
tonian simulation. The increase of the energy content in 
the flow is even more apparent if we consider that, in 
addition to kinetic energy, a non-negligible contribution 
arises from the free energy a. In terms of total energy, 
E = 3/2 Ur ms +a, the increase with respect to Newtonian 
turbulence is of the order of 48%. Clearly, the increase 
in the total dissipation for fixed total viscosity amounts 
to a reduction of the conventional Kolmogorov scale. 

A further effect is observed when comparing the Tay- 
lor microscale for the two simulations. The increase in 
turbulent kinetic energy of about 44% is accompanied by 
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200 400 
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FIG. 2. Top: Time behavior of the integral scale for the 
Newtonian, (dotted) and viscoelastic, (solid) case. Bottom: 
Time behavior of the spatial average of the turbulent kinetic 
energy, same symbols. 
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FIG. 3. Top: Time behavior of the spatial average of 
the enstrophy for the Newtonian, (dotted), and the viscoelas- 
tic case, (solid). Bottom: Time behavior of the Taylor mi- 
croscale, same symbols. 

a reduction of the order of 50% in the enstrophy, giving 
an increase of 73% for the Taylor microscale. This sug- 
gests a substantial depiction in the energy content of the 
high wave- number modes. Consistently a substantial in- 
crease in the Taylor-Reynolds number, Re\ = u rms \jvT-, 
is achieved. 

From the numerical point of view, a comparable value 
for Rex in Newtonian fluids cannot be reached with only 
64 spatial modes. For the viscoelastic case, however, the 
requirements on the grid are less strict, since most of the 
dissipation, the polymeric contribution, is not associated 
to spatial gradients. We observe in passing, that spectral 
calculation with the FENE-P model may be prone to 
high-wavenumber instability, which is maintained under 
control by adding a small amount of artificial viscosity in 
the algorithm for the conformation tensor see also 
]l9[ for an alternative approach. 

The evolution of the system, followed for about 
50 large-eddy turnover times, T = (Ix/fo) 1 ^ 2 = 
7.9, is globally described in figure ||, where the 
time history of the instantaneous integral scale L — 
ir/2 J Q mam l/ke(k) dk / v% ms , where e is the instanta- 
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neous energy spectrum, is shown on the top. From the 
comparison between viscoelastic and Newtonian case, the 
average increase in L is apparent. Clearly this is associ- 
ated to the increase in the amount of energy in the energy 
containing range, and corresponds to the observation of 
larger structures in the field. The presence of fewer but 
more energetic coherent structures is also consistent with 
the existence of more pronounced intermittency cycles, 
more evident in the history of the kinetic energy shown 
on the bottom part of the same hgure. 

Figure [| describes the evolution of enstrophy, top, 
and Taylor microscale, bottom. On average, the enstro- 
phy reduces with respect to Newtonian turbulence, and 
presents large spikes following with slight delay the cor- 
responding peaks found in the kinetic energy. In figure ^ 
we report on top the time behavior of the spatial aver- 
age of the free energy. Apparently, the amount of energy 
stored by the polymers is substantially lower than the 
kinetic energy, already given in figure |2| and reproduced 



X, 1.5 





400 



FIG. 4. Viscoelastic simulation, De — .38. Top: Time 
behavior of the spatial average of the free-energy, a, (solid) 
compared with the kinetic energy (dotted). Note the scales, 
differing by an order of magnitude, for free-energy on the right 
and kinetic energy on the left. Bottom: Time behavior of the 
spatial average of the Newtonian (dotted) and the viscoelastic 
component (solid) of the dissipation. 



here for convenience. From the comparison one may 
note the more irregular behavior of the free-energy with 
reduced peak to peak oscillations. Finally, on the bottom 
of the figure we plot the evolution of the two forms of dis- 
sipation, Newtonian and viscoelastic respectively. Here 
the trend is opposite, with a larger portion of dissipation 
occurring in the polymers. 



VI. SCALE BY SCALE BUDGET 

The Karman-Howarth equation ( |l6| ) provides the ap- 
propriate tool to analyze the scale by scale budget for the 
velocity field. To ease its interpretation one should keep 
in mind the term 



1 d 

nr 2 dt J B 



Ca dV r 



originally present at the left-hand side of the equation, 
that has been dropped due to the assumed stationarity 
of the field. 

A negative flux through the boundary of B r , 



$(r) 



f 



Airr 2 



dB r 



< SV 2 5V k + 2 T fc * 5Vi > n k dS r , (24) 



corresponding to a flux towards the inside of dB r , implies 
a negative contribution to the integral up to r of the 
correlation, i.e. an average depletion of correlation in the 
scales smaller than r. For stationary flows, the depletion 
balances the feeding operated by the forcing, 



F(r) 



-4/3 < f iUi > r = -4/3e T r, 



consistently with a steady correlation at all scales. 

For polymers, the flux is split into two parts, the classi- 
cal part which, for isotropic conditions, can be expressed 
in terms of the third order longitudinal structure func- 
tion, 



$ c (r) 



47rr 2 



< SV 2 SV k > n k dS r = - < 5V? > 



0B r 



(25) 



and the additional term associated with the viscoelastic 
reaction, 



1 



Anr 2 



<T* l 5V l >n k dS r 



(26) 



dB r 



Figure [| shows all the contributions appearing in the 
steady state Karman-Howarth equation, and, in partic- 
ular, on the top, it reports the budget of the Newto- 
nian simulation for comparison. Considering the Newto- 
nian case first, the straight-line with negative slope corre- 
sponds to — 4/3eArr, filled squares. The positive correla- 
tion between forcing increments and velocity increments 
is represented by the filled diamonds. The nonlinear 
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transfer term, <I> C , is given by the filled circles, while the 
viscous correction corresponds to the filled triangles. All 
terms sum up to zero within the accuracy of the available 
statistics. Clearly, the viscous correction is present at all 
scales, confirming that the Reynolds number of the sim- 
ulation is too small to achieve a proper inertial range. A 
further point to be highlighted is the non-negligible con- 
tribution of the correlation between forcing increments 
and velocity increments. Nonetheless, for a fair range 
of scales, up to r/rjx ~ 30, the balance involves only 
the nonlinear transfer term , <E> C , and the forcing F. $ c 
displays a behavior in terms of separation which, a part 
from minor corrections due to the other terms, follows 
the linear trend of F. 

The same kind of budget is discussed for the viscoelas- 
tic case in the bottom part of the figure. The filled sym- 
bols have the same meaning defined for the Newtonian 
case, except that, now, the transfer term, filled circles, 
corresponds to the entire flux i.e. it accounts also for 
the polymeric contribution. As we may see, the picture 
grossly reproduces that already described for Newtonian 
turbulence. One should note that the forcing term, F, 
now involves the entire dissipation ex, i.e. the sum of the 
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FIG. 5. Scale by scale budget for the Newtonian (top) 
and the viscoelastic case (bottom). Symbols are defined in 
the text. 



Newtonian and the polymeric contributions. The slope 
is then substantially larger than in the Newtonian case 
displayed on the top. Not unexpectedly, the global flux 
term is qualitatively similar to the corresponding New- 
tonian one, with a substantial range of scales where the 
linear trend defined by F is reasonably well reproduced 
by $. This emerges quite clearly from the plots in fig- 
ure [j] where the budgets for Newtonian and viscoelastic 
turbulence are superimposed the one to the other, after 
a rescaling proportional to the respective values of the 
global dissipation. Also for the polymers, we are in pres- 
ence of a direct cascade, with a negative global flux. 

Returning to the plots on the bottom of figure [| the 
most interesting feature is represented by the open sym- 
bols, circles and squares, which give the splitting of the 
flux into the classical contribution associated to the third 
order velocity structure function, $ c (circles), and to the 
polymer contribution, <& p (triangles), respectively. The 
larger scales are less affected by the polymers, with $ p 
becoming sub- leading with respect to «3> c . One may con- 
jecture that, at sufficiently large scale and for Reynolds 
numbers large enough, the turbulence scalings are sub- 
stantially unaffected by the polymers. Under this as- 
sumption, we see that the total dissipation ct may be 
taken to measure a purely inertial energy transfer at large 
scales. This attaches to our conventional definition of 
Kolmogorov scale the obvious meaning of the length-scale 
where dissipation would occur in the absence of polymers. 

Within the scenario according to which the large scales 
are unaffected by the polymers, the time criterion of 
Lumley [Q fixes the scale r* below which the poly- 
mers begin to feel the turbulent fluctuations. This hap- 
pens when the eddy-turnover time of classical turbulent 
eddies, i.e. without the effect of the viscoelastic reac- 
tion, equals the principal relaxation time of the polymer 
chains. Assuming Kolmogorov scaling for the velocity 
increments, this leads to 

r* = (e T r 3 ) 1 / 2 (27) 

or, scaled with the conventional Kolmogorov length, to 

r*/vr = r 3 / 2 (e T /v T ) 3/4 ■ (28) 

From the data given in table Q, for our simulation we have 
r* /ifr — 70, within a computational box of l x /r]T — 174, 
suggesting that the forcing, applied to the first shell of 
wavenumbers, acts just above the scale where the poly- 
mers already feel the turbulent fluctuations. In fact, a 
better way to apply the time criterion is to identify the 
eddy-turnover time as function of separation from the ac- 
tual data provided by the DNS. This is done in figure [?], 
where as characteristic time for the scale r we assume 

1 /3 

r/<P c <. From the fi gure we infer that our viscoelastic tur- 
bulence achieve an eddy-turnover time comparable with 
the viscoelastic relaxation time, horizontal line in the fig- 
ure, at a scale considerably smaller then predicted on 
purely dimensional grounds, leaving a reasonable range 
of between the forcing and Lumley scale. 
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Coming back to our comparison between the two com- 
ponents of the flux, at small scales we find that the com- 
ponent originated by the polymers takes over, entailing 
the reduction of the third-order structure function. The 
two curves, open circles and squares, clearly identify a 
cross-over scale £ p according to the equation 

$ P (£ P ) = $ c (£ p ) . (29) 

For our case the balance is achieved at r/r/T — 22, while 
\/rjT = 33. Observe for comparison that, in the corre- 
sponding Newtonian case, A/77 is about half the value. 
By inspection of the figure, we also find that all along 
the available range of scales the polymers reaction is dy- 
namically relevant for the turbulent fluctuations, with $ p 
always comparable in magnitude with <E> C . 
Dimensionally, the right-hand side of eq. (^9|) can be 
estimated as It £ P , as follows by roughly assuming a 
purely classical inertial scaling for r > £ p . On the 
other hand, $ p (£ p ) should be somehow proportional to 

a** v p jr e^ 3 £ P ^ 3 , where, again, Kolmogorov scaling is 
used to estimate 8V in eq. (|2q). Accounting for the con- 
stitutive relation (Q), T* is expressed in terms of the, as 
yet, undetermined factor a** which should account for 
the level of stretching achieved in the polymers at the 
cross-over scale. This yields the expression 

*p _ 1 

which can be used to build the dimensionless quantity 

a** = rjv p (£ p e T ) 2 / 3 . (30) 

For the present case we have (v p /t) 3 ' 2 /€t — 1.4 10~ 5 , 
which, after considering that £ p /i] p ~ 22, 



yields a** ~ 1400. The parameter a** allows to quan- 
tify the effective stretching of the chains at cross-over by 
assuming a** oc f**R* k * k /pl, see eq. (§). 

The quantity 

& = f{Rkk', Pmax, Po)Rkk/ Po (31) 

is plotted versus the elongation Rkk/po on the right of 
figure |^. We see clearly the finite extensibility effect, 
with the unbound growth of a as the critical extension 
Pmax/Po i s approached from the left. Note that, at small 
elongation, a goes almost linearly with R k .k/Po- The hor- 
izontal line plotted in the figure corresponds to the value 
a** we have estimated from the numerical simulation. 
The intercept of this line with the curve representing a 
lays relatively far from the critical extension, though well 
within the region of non-linear elastic response. 

It may be instructive to compare the elongation asso- 
ciated to the cross-over scale to the pdf of Rkkj p\- 
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FIG. 8. a, see eq. ([u]), vs. polymers elongation, Rkk/po- 

To this purpose, we superimpose to the plot of a given 
in figure || the histogram of the elongation as obtained 
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by the numerical simulation. The effective elongation at 
cross-over, denoted by the vertical line, lies to the right 
of the maximum of the histogram, indicating a level of 
stretching substantially larger then either the most prob- 
able or the average elongation. Clearly the estimate here 
proposed should be understood as an order of magnitude 
analysis, to allow the comparison of different flow condi- 
tions. 

To complete the analysis, in figure ^| we address the 
Yaglom-like equation for the free-energy. As seen in the 
top part of the figure, the left and the right hand side 
of eq. (|2^) balance within the accuracy of the available 
statistics. 

The fluctuations of excess-power are converted into a 
free-energy flux with a negative sign, as wee see from 
the figure. Hence, also for the microstructure, we are in 
presence of a direct cascade. We like to stress that the 
direction of the cascade is not a priori obvious in this 
case, since we miss the standard term, linear in the 
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FIG. 9. Top: The free-energy flux (filled triangles), and the 
forcing term given by the excess-power (filled circles) accord- 
ing to Yaglom equation. Bottom: The two contributions to 
the excess-power term in the Yaglom equation, stress-power 
term (filled circles) and polymeric dissipation term (filled tri- 
angles) . 



separation and related to the dissipation of the scalar, 
which usually determines the sign of the energy flux. 
From the figure, the free-energy flux, 



< Sa 2 SV\\ > . 



manifests a pronounced minimum in correspondence with 
the cross-over scale £ p , suggesting that, near cross-over 
the activity of the polymers is particularly intense and 
the polymers tend to become less active as the scale is 
increased. The bottom part of the figure gives the de- 
composition of the right-hand side of equation (|22|), with 
the positive term 

associated to the stress power denoted by filled squares 
and that associated with the polymeric dissipation indi- 
cated by the filled triangles. 

Order of magnitude considerations confirm that the 
data shown in the bottom part of figures ^| and in figure 9 
are consistent. At cross-over, r = £ p , we find from figure 5 
that $** ~< SV** 3 >~ .12, suggesting SV** ~ .5. From 
figure § we also have < Sa** 2 SV** >~ .00015. It follows 
Sa** ~ .017 a value considerably close to SV** 2 ~ .025. 
In other words, at the cross-over scale fluctuations of 
free-energy, or more loosely of elastic energy, and fluc- 
tuation of turbulent kinetic energy become comparable, 
as it should be expected. 



VII. CONCLUDING REMARKS 

The interaction between long chain polymers in dilute 
solution and turbulence is known to have a profound in- 
fluence on the structure of turbulence. Several possible 
explanations have been proposed, with the eventual aim 
to understand the drag reducing properties manifested 
by such fluids in wall bounded turbulent flows. Even 
the most powerful approach, however, has to face the 
elusive mechanism by which polymers dynamics and tur- 
bulence are mutually coupled. Among others, specially 
noteworthy is the theory of De Gennes and Tabor, who 
proposed a cascade model for the viscoelastic turbulence 
of dilute polymers solutions. Their elegant theory is- 
sues from considering homogeneous isotropic conditions 
and rests upon suitable assumptions on the self-similar 
stretching of polymers chains under a fluctuating velocity 
field. The entire model is essentially phenomenological in 
nature, and should thus be contrasted with experimental 
observations. 

In this context, we believe to have provided here a 
reasonable framework based on two essential ingredients. 
The first ingredient is the direct numerical simulation 
of the natural flow where the cascade issue should be 
addressed, namely homogeneous isotropic turbulence nu- 
merically modeled as a triply-periodic box with random 
forcing at large scales. Beside the obvious implementa- 
tion in terms of spectral algorithms, this step required 
the selection of a good rheological model for the poly- 
mers. Based on previous experience, our natural choice 
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has been the FENE-P model, known to capture reason- 
ably well the drag-reducing effects still maintaining an 
affordable level of computational complexity. 

The second ingredient is the extension to viscoelastic 
fluids of two of the, let's say, few exact results known in 
classical turbulence theory, namely the Karman-Howarth 
equation for the velocity and the Yaglom equation for 
scalars. As a characteristic feature of viscoelastic tur- 
bulence, the coupling between velocity and extra-stress 
in the Karman-Howarth equation prevents recasting the 
result in the familiar form involving only longitudinal 
structure functions. Hence we have adopted a natural 
generalization in terms of suitable surface and volume 
integrals in the space of separations. Concerning the mi- 
crostructure, the descriptor field is a second order ten- 
sor, making the derivation of the appropriate equation 
exceedingly cumbersome. We decided then to fucus our 
attention on the single most important scalar quantity, 
namely the free-energy, whose Yaglom-like equation fol- 
lows straightforwardly. 

The scale by scale budget based on the Karman- 
Howarth and Yaglom equation has been used to under- 
stand the respective role of the traditional energy transfer 
term versus polymeric transfer. We find a direct cascade 
occurring both in the kinetic field and in the microstruc- 
ture, for which, as already discussed, the sign of the en- 
ergy flux is not a priori obvious. 

The budget has allowed us to identify a cross-over scale 
below which the polymers contribution becomes domi- 
nant, being sub- leading with respect to incrtial transfer 
at larger scales. This is the analogous of the scale defining 
the elastic limit in the theory of De Gennes and Tabor, ft 
is here derived in a rigorous context, and evaluated on the 
basis of a direct numerical simulation of a drag-reducing 
visco-elastic material as described by the FENE-P model. 

We find that the polymers substantially deplete the 
energy content of the smallest scales, in agreement with 
recent experimental results, but they affect the dynam- 
ics of the fluctuations at all the scales of our simulation. 
The latter effect could be interpreted as an artifact of 
the limited extension of our computational domain which 
forced us to apply the external excitation relatively close 
to the scale defined by Lumley's time criterion. Nonethe- 
less the present results make questionable the existence 
of a purely passive range, where polymers are deformed 
without back-reaction on the velocity field. 

Based on the balance of terms appearing in the exact 
form of the Karman-Howarth equation we have proposed 
a dimensionless parameter which may give a quantita- 
tive measure of the stretching the polymers experience at 
cross-over, when elastic effects become comparable with 
the corresponding kinetic energy. 

Finally, a substantial level of intermittency is observed 
in the system, with fluctuations in turbulent kinetic en- 
ergy considerably larger than usually observed in homo- 
geneous isotropic turbulence of Newtonian fluids. This 
effect is presumably associated with the existence of the 
additional mechanism of energy removal from the scales 



of the forcing. In particular, a substantial part of the 
energy income does not follow the classical route cascad- 
ing towards viscous dissipation. Instead, it is moved to 
the microstructure, to feed an additional cascade process 
seemingly requiring a larger level of fluctuations. 
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